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Abstract 

The shape and tails of partial distribution functions (PDF) for a climatological signal, i.e. the El 
Nino SOI and the turbulent nature of the ocean-atmosphere variability are linked through a model 
encompassing Tsallis nonextensive statistics and leading to evolution equations of the Langevin 
and Fokker-Planck type. A model originally proposed to describe the intermittent behavior of 
turbulent flows describes the behavior of the normalized variability for such a climatological index, 
for small and large time windows, both for small and large variability. This normalized variabil- 
ity distributions can be sufficiently well fitted with a ^^-distribution. The transition between the 
small time scale model of nonextensive, intermittent process and the large scale Gaussian exten- 
sive homogeneous fluctuation picture is found to occur at above ca. a 48 months time lag. The 
intermittency exponent (k) in the framework of the Kolmogorov log-normal model is found to be 
related to the scaling exponent of the PDF moments. The value of /t (= 0.25) is in agreement with 
the intermittency exponent recently obtained for other atmospheric data. 
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I. INTRODUCTION 



Fractional Gaussian noises and fractional Brownian motions [1] have served recently as 
models for a wide variety of data in various fields, like meteorology [2, 3], geology [4], cardiac 
dynamics [5], finance [6-8]. The concept of fractional Gaussian noise (fGn) as formulated 
by Mandelbrot and Van Ness [1] is a derivative process obtained from fractional Brownian 

motion (fBm) Buit)^ namely lim^ .i^Bnit + 5) — BH{t))/5. As it has been shown by 

Flandrin [9] . Even though the fBm is a non-stationary process it obeys a power law over all 
frequencies. In this study we use fGn as a model for a climatological signal. 

The origin of non-Gaussian, thus non negligible large and sometimes so called extreme 
volatility events characterized by so called fat tailed distributions is a key question in sta- 
tistical physics; the fat tails of (short and long-range) volatilities are thought to be caused 
by some 'dynamical process'. Destroying all correlations, e.g. by shuffling the order of the 
fluctuations, is known to cause the fat tails almost to vanish. 

The fat tails indicate an unexpected high probabihty of large changes. These extreme 
events are of utmost importance for risk analysis. They are considered to be a set of strong 
bursts in the energy dissipation. In so doing the PDF and the fat tail event existence are 
thought to be similar to the notion of intermittency in turbulent flows [10]. 

It is an open question whether both the fat-tailed power-law of partial distribution func- 
tions (PDF) of the various volatilities and their evolution for different time delays in clima- 
tological indices can be described. 

On the other hand, the non-Gaussian character of the fully developed turbulence [11] has 
been linked to nonextensive statistical physics [12-20]. It seems that there is no study of 
TsaUis statistics application or approach in chmatology [21]. 

One of the most intriguing phenomena in climatology, known as El Nino, i.e. the more 
or less cyclic warming and cooling of the eastern and central regions in the Pacific Ocean, 
appears to be very complicated to describe [22-25]. There are three sorts of models, based 
on special types of filters designed specifically to detect a signal from given atmospheric 
input. The abihty of these models for producing successful forecasts of El Nino appears 
to be concomitant to the very low-frequency and large-scale evolution of the characteristic 
patterns in the atmospheric boundary layer circulation. This evolution can be thought of as 
a signal that precedes El Nino events. It would be useful to have model-independent results 
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with basic geophysical inputs. 

El Nino is a disruption of the ocean-atmospheric system in the tropical Pacific having 
important consequences for weather around the globe. It is factually described by the so 
called Southern oscillation index (SOI). Much of the drastic and tragic events occurring 
in North America, Tropical Africa and Europe are also attributed to El Nino. There are 
different ways of characterizing El Nino events. One of them is by the values of an index 
such as the Southern Oscillation Index (a proxy measure of El Nino based on surface air 
pressure differences between Darwin, Australia and Tahiti, French Polynesia) [29], or large 
sea surface temperatures (SST) and sea surface height (SSH) anomalies in the eastern trop- 
ical Pacific Ocean. Here we will use the sea level pressure (SLP) differences between two 
meteorological stations, one at Tahiti, the other at Darwin. Questions on variability of the 
SOI on various time scales are relevant for better modelling. It was recently found for the 
southern oscillation index (SOI), characterizing El Nino events [26] that long-range correla- 
tions exist between the fiuctuations of the index. Also, correlation between the SOI and the 
North Atlantic Oscillation (NAO) index has been reported [27] having a common oscillation 
of about 6 — 8 years. The NAO was studied per se, in [28]. 

In this paper the behavior of a climatological index, i.e. the SOI, on short and large 
time windows (time scales to be better defined below) is studied along the lines of a re- 
cently suggested model of hydrodynamic turbulence that serves as a dynamic foundation for 
nonextensive statistics [16-18]. 

In Sect. 2, we describe the distribution of variability for the monthly value signal of the 
SOI index for the time interval between Jan. 1866 and Jan. 2006, thus a series of N =1681 
months or data points, values downloaded from the Climate Prediction Center web site 
{http : / /www. cpc.ncep.noaa.gov /data/ indices/). We also characterize the tail(s) of the 
distribution for various time lags Ai's, i.e. from 1 month up to 36 months (or 3 years!), 
and will observe the value of the PDF tails, for such time lags, outside the best gaussian 
fit through the data. There are several ways of displaying features in the variability of a 
climatological index. A simple one represents value increment Ay{t) = y{t + At) — y(t) or 
difference between the value of the index y{t) at time t + At and its value at time t. Below 
we mainly consider the normalized variability Z{t,At) = {Ay{t)— < Ay >At)/o"Ai- where 
< Ay >At denotes the average and (jAt the standard deviation of Ay{t) for a given At. 
The normahzed variabihty Z{t, At) depend on the time t and the time lag At. However, 
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in order to simplify the notations and whenever possible without leading to confusion and 
misunderstanding we will drop the explicit writing of one or both variables. 

In Sect. 3, we calculate the power law exponents characterizing the integrated distribu- 
tion of the normalized variability over different time lags for the SOI index monthly values 
through a detrended fluctuation analysis and a power spectral density analysis point of view. 

In all cases it is useful to test the null hypothesis or estimate the error bars with respect 
to standard signals. It is thought [10] that the fat tails are caused by long-range volatihty 
correlations. Destroying all correlations by shuffling the order of the fluctuations, is known 
to cause the fat tails almost to vanish. A Kolmogorov-Smirnov test (not shown) on shuffled 
data has indicated us the statistical validity of the numerical values and the statistically 
acceptable meaning of the displayed error bars. Results are compared to shuffled data for 
estimating the value of the error bars. 

In Sect. 4, Tsallis statistical approach is outhned, and distributions of (normahzed) 
variability for time lags between At = 1 to 36 months are examined. The nonextensivity, 
i.e. some anomalous scaling of classically extensive properties like the entropy, is linked to 
a single parameter q, e.g. in the Tsallis formulation of nonextensive thermostatistics. It 
is found that the g- value of the nonextensive entropy converges to a value = 1.01 for At 
— 36 months, starting with q — 1.25 for At — 1. The probability density fAt{P) of the 
volatility (3 in terms of the standard deviation of the normalized variability of the SOI for 
different time lags is found to obey the x^-distribution. The intermittency exponent (k) of 
the Kolmogorov log-normal model is found to be related to the scahng exponent of the PDF 
moments, -thereby giving weight to this model. 

In Sect. 5, the usual Fokker-Planck approach for treating the time-dependent proba- 
bility distribution functions is summarized. Coefficients governing both the Fokker-Planck 
equation for the distribution function of normalized variability and the Langevin equation 
for the time evolution of normalized variability of monthly value signal of SOI are obtained. 
Therefore we present for the first time a coherent theory linking the shape and tails of partial 
distribution functions for long and short time lags of the monthly values of a climatological 
signal and connect the often suggested turbulent nature of the ocean-atmospheric interface 
interactions to a model encompassing nonextensive statistics and evolution equations of the 
Langevin and Fokker-Planck type. 

We will often compare results based on normalized variability and non normalized vari- 
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ability time series. 



II. DATA AND DISTRIBUTION OF VARIABILITY 

Monthly values of the SOI index for the time interval between Jan. 1866 and Jan. 2006 
were downloaded from the Joint Institute for Study of the Atmosphere and Ocean (JISAO) 
web site http : // tao.atmos.washington.edu/ pacs/ additional ^analyses/soi. html for the 
longest period available there, i.e. from Jan. 1866 to June 1999. Data for the time interval 
from July 1999 to Jan. 2006 were downloaded from the Climate Prediction Center NCEP 
web site {http : / /www. cpc.ncep.noaa.gov /data/ indices/). For the years before 1866 daily 
measurements of the sea level pressure at both stations have been reported to exist and 
monthly values of the southern oscillation index have been calculated [34] back to 1841. 
However, there are gaps of a couple of years in the record whence not suitable for our 
analysis. We have chosen the JISAO's data for the period before July 1999, since they do 
not contain missing values as the data series from NCEP do. The data are fully compatible 
because the standardization of both JISAO and NCEP data series is calculated through the 
standard deviation S of the sea level pressure (SLP) at a station in Tahiti and the sea level 
pressure at a station in Darwin 

p p 

Tahiti -'■ Darwin ^-^\^ 

S 

Such SOI monthly data are plotted as a function of time in Fig. 1. The data set consists 
of 1681 data points. It is sometimes stated that daily or weekly values of the SOI do not 
convey much in the way of useful information about the current state of the climate, and 
accordingly the Bureau of Meteorology does not issue them. Daily values in particular can 
fluctuate markedly because of daily weather patterns, and should not be used for climate 
purposes. We may disagree with this statement (provided the reliability of the data). There 
are indeed techniques which can sort out noise from coherent behavior [30] . 

Sustained negative values of the SOI often indicate El Nino episodes [31]. These negative 
values arc usually accompanied by sustained warming of the central and eastern tropical 
Pacific Ocean, a decrease in the strength of the Pacific Trade Winds, and a reduction in 
rainfall over eastern and northern Australia. The most recent strong El Nino was in 1997/98. 
Positive values of the SOI are associated with stronger Pacific trade winds and warmer sea 
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temperatures to the north of Austraha, popularly known as a La Nina episode. Waters in 
the central and eastern tropical Pacific Ocean become cooler during this time. Together 
these give an increased probability that eastern and northern Australia will be wetter than 
normal. The most recent strong La Nina was in 1988/89; a moderate La Nina event occurred 
in 1998/99, which weakened back to neutral conditions before reforming for a shorter period 
in 1999/2000. This last event finished in Autumn 2000. 

The distribution of the normalized variabihty Z{t, At) of the monthly value signal of SOI 
index between Jan. 1866 and Jan. 2006, for At — 1 month are plotted in Fig. 2a. The 
partial distribution of non normalized variability Ay of the monthly value signal of SOI 
index is plotted in Fig 2b as we will later compare some of our findings for both Z and Ay. 
For comparison A fit is first attempted with a Gaussian distribution for small values of the 
increments, i.e. the central part of the distribution. The distribution is well fitted with such 
a Gaussian type curve within the interval Z e] — 2, 2[ but departs from the Gaussian form 
outside this interval. The negative and positive tails of the distribution outside the Gaussian 
curve are found both to be equal to -5.6. In the case At > 1 month, it is observed that the 
best Gaussian range is increasing with increasing time lag (Fig. 3). 

III. TIME CORRELATIONS AND SPECTRAL POWER 

There are different estimators for the long and/or short range dependence of fiuctuations 
correlations [36]. 

Through the (linearly) detrended fluctuation analysis (DFA) method, see e.g. [37], we 
show first that the long range correlations of monthly value signal of SOI for the time 
interval of interest, are l//-like. The method has been used previously to identify whether 
long range correlations exist in non-stationary signals, in many research fields such as e.g. 
finance [7, 8], cardiac dynamics [5] and of course meteorology [2, 3, 38]. Its concepts are 
therefore not repeated here. For an extensive fist of references see [37]. Briefly, the signal 
time series y{t) is first integrated, to 'mimic' a random walk Y{t). The time axis (form 
1 to A^) is next divided into non-overlapping boxes of equal size n; one looks thereafter for 
the best (linear) trend, Zn, in each box, and calculates the root mean square deviation of 
the (integrated) signal with respect to Zn in each box. The average of such values is taken 
at fixed box size n in order to obtain 
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F{n) 



1 ^ 

-J2[Y{t)-z4t)f (2) 

i=l 

The box size is next varied over the n value. The resulting function is expected to behave 
like n^^^^-^-^ indicating a scaling law. For the (integrated) monthly value signal of the SOI 
index, a scaling exponent 1 + Hofa = 1-05 ± 0.01 is found (Fig. 4) in a scaling range 
extending from about 4 to about 66-72 months. A signal with Hausdorff dimension Hdfa 
close to zero has the characteristics of a fractional Gaussian noise signal [4, 39]. 

Along the same line of thought the scaling properties of the normalized variability 
Z{t,At) = {Ay{t)— < Ay >At)/o'At have also been tested for different time lag values, 
i.e. At = 1,3,6,12,24,36 months (Fig. 5). The DFA functions, as defined here above, of 
the integrated normalized variability shows non trivial scaling properties for the series of 
normalized variability. The values of the scaling exponents and the maximum box size 
(in days) for which the scaling holds for each DFA-function are given in Table I, while the 
DFA-functions together with fitting lines are plotted in Fig. 5. 

The power spectrum of the monthly value signal of SOI S{f) ~ with spectral expo- 
nents fii = —0.26 and ^2 = 1-20 with a scale break at 1/70 months^^ is shown in Fig. 6. The 
scaling properties of the power spectrum of two surrogate data, one in which the amplitudes 
are randomly shuffied and another in which the magnitudes are preserved but the sign of 
the data is shuffled, are shown in the inserts of Fig. 6. Such scahng spectral exponents 
II — are signature of a white noise like behavior. Recall that /i — 2.0 corresponds to usual 
Brownian motion. The theoretical relationship /i — 2Hdfa — 1 is approximately verified, - 
the weak agreement being likely due to the quite limited data size. 

We have also checked for scaling behavior and possible periodicities in the power spectrum 
of the time series of the normalized variability Z{t,At) — {Ay{t)— < Ay >At)/cAt for 
different (selected) values of the time lag At — 1, 3, 6, 12, 24, 36 months (Fig. 7). 

Periodicities in the power spectrum of the normalized variability time series for At > 
1 month were expected to be found since these periods are somewhat embedded into the 
time series by the way they are obtained and the Fourier transform technique. It is easily 
observed that the maxima and the minima of the spectrum correspond to harmonics and 
subharmonics of 1/At. 
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IV. TSALLIS STATISTICS 



Based on the scaling properties of multifractals [40] Tsallis [12, 41] proposed a generalized 
Boltzmann-Gibbs thermo-statistics through the introduction of a family of non-extensive 
entropy functional Sg given by: 

Sg = ( 1 - / p{x, tydx] , (3) 



q-1 

with a single parameter q and where /c is a normalization constant. The main ingredient in 
Eq.(3) is the time-dependent probability distribution t) of the stochastic variable x. The 

functional is reduced to the classical extensive Boltzmann-Gibbs form in the limit of q — > 1. 
The Tsallis parameter q characterizes the non-extensivity of the entropy. Subject to certain 
constraints the functional in Eq.(3) seems to yield a probability distribution function of the 
form [10, 12, 16, 32, 33] 

,w ' /i^g/^f°(^-'yr'^ (4) 



Zg[ 2a -{q- 1) 

for the stochastic variable x, where 



f CPo2a{q - 1) ) r (^) 



2Q;-(g-l) j Y(^)r(^-^ 



in which C is a constant and < a < 1 is the power law exponent of the potential U {x) = 
C\x\'^" that provides the 'restoring force' F{x) in Beck model of turbulence [16-18, 20]. The 
latter is described by a Langevin equation 

^ = -7F(x) + R{t) (6) 

where 7 is a parameter and R(t) is a gaussian white noise. A non-zero value of 7 corresponds 
to providing energy to (or draining from) the system by the outside [42] . The parameter /3o 
in Eq.(4) and (5) is the mean of the fluctuating standard deviation /3, i.e. the local standard 
deviation of \x\ over a certain window of size m [10]. We will use this model assuming that 
the normahzed variabihty Z{t, At) represent the stochastic variable x, as in Eq.(l). We will 
search whether Eq.(4) is obeyed for x = Z{t,At), thus studying p{x) = pAt{Z) for various 
time lags At. 
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Just as in Beck model of turbulence [16-18] we assume that the standard deviation (5 is 
X^-distributed with degree v (see another formula in [20]): 

/a*(/5)- w^f;^y^%^/^-^expf-^V ^ > 2, (7) 



T{v/2) \2PoJ " "V 2/3o 

where F is the Gamma function, /3o =< P > and the number of degrees of freedom u can be 
found from: 

2 < /? >2 



< /52 > - < 

The TsaUis parameter q satisfies [16] 



(8) 



2q; 

q = l + ——. 9 
aiy + 1 

To justify our assumption that the 'local' standard deviation of the normalized variability 
Z{t, At) is of the form of x^-distribution, we checked the distribution of the normalized 
variability of the monthly value signal of SOI. We have calculated the standard deviation of 
the normalized variability within various non-overlapping windows of size m, ranging from 
6 to 36 months 



m 



(fe+l)m / (fe+l)m \ 

i=km+l \ i=km+l J 

In doing so we have a various number of M. non-overlapping windows for various time 
lags At, and have searched for the most efficient size of the window in order not to loose 
data points and therefore, information. The resulting empirically obtained distributions of 
the 'local' standard deviation (Eq.(lO)) of normalized variability for the difi^erent time lags 
of interest are plotted in Fig. 8 for an intermediary case m = 12. The values of the degree 
u of the x^.fjistribution are then obtained using Eq. (8). The spread [Pmin, Pmax] of the 
local standard deviation /3 decreases with increasing the time lag as it is expected from a 
X^-distribution function due to the exponential function in Eq. (7) for large values of the 
degree of freedom i/. The value of u much varies as a function of m and the time lags 
considered. The fits arc always very good. However the and u values are quite dependent 
on the parameters used in the numerical analysis. Based on these results, e.g. Fig. 8, it can 
be accepted that the (turbulent) model /^-distributions can be sufficiently well fitted for our 
purpose with a x^-distribution, thereby justifying the initial assumption. [57] 
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The probability distributions of the normahzed variabihty for the different values of the 
time lag At — 1, 3, 6, 12, 24, 36 months are shown in Fig. 3 together with the lines repre- 
senting the best fit to the TsaUis type of distribution function. In Table II the statistical 
parameters related to the Tsallis type of distribution function are summarized, including a 
criterion for the goodness of the fit, i.e. the Kolmogorov-Smirnov distance d^s^ which is de- 
fined as the maximum distance between the cumulative probability distributions of the data 
and the fitting hues. Note that the kurtosis (see Table II) for the TsaUis type of distribution 
function 

(U) 

where Kl = 3 for a Gaussian process, is positive for all values of g < 7/5 as expected, since 
its positiveness is directly related to the occurrence of intermittency [10]. Moreover, the 
limit q < 7/5 also implies that the second moment of the Tsallis type distribution function 
will always remain finite, as necessarily due to the type of phenomena hereby studied. 

Furthermore, if we assume that the Kolmogorov log-normal model of turbulence [47] is 
applicable and let Ati be the scale at which the whole partial distribution function becomes 
Gaussian, then the kurtosis Kr should scale as 

if, = X.(|i)". (12) 

Therefore 

5-7(At/At,r (13) 

3-5(At/Ati) 

In order to obtain an estimate for At^, wc increase the time lag to the value At = 
48 months, quite outside the range so far examined (see Fig. 3 for example) leading to a 
rather complete coincidence between the distribution functions in the Tsallis and Gaussian 
forms for the presently investigated data. The corresponding parameter values are listed 
in Table II. A quick perusal observation convincingly indicates where the transition occurs 
between the small time scale model of nonextensive, intermittent process and the large scale 
Gaussian extensive homogeneous fluctuation picture [10, 12]. 

In Fig. 9 the Tsallis parameter q is shown as a function of the rescaled time lags At/ At^, 
where Ati is the integral scale, the scale at which the whole probability distribution function 
converges to Gaussian. The crosses represent the q values for which the best fit to the SOI 
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data (Fig. 3) is obtained with Eq. (4). With this the value of the integral scale Ai^, we find 
the value of the exponent S — 0.11 as the one for which the Eq. (13) fits best the g'-values. 
The exponent value 5 = 0.11 also allows to fit well the power law dependence (Eqs. (11) 
and (12)) of the rescalcd kurtosis K,,./ as shown in the insert of Fig. 9. 

Note that in the framework of the Kolmogorov log-normal model [17, 47], 5 = 4k/9, where 
K is called the intermittency exponent. Therefore, we find k — 0.25 for the intermittency 
exponent of normalized variability of the SOI signal in the time interval of interest. This 
value of K is interestingly the same as the value of the intermittency exponent k — 0.25 
for turbulence recently obtained from experimental atmospheric data [48]. Early estimates 
have varied from 0.18 to 0.85 using different experimental techniques [49-51]. Large range of 
values of the intermittency exponent, ranging from 0.2 to 0.8, have been reported in studies of 
multiparticle production [52] . It was found that the range of intermittency exponent values 
depend on the number of cascades; the smaller the number of stages of the multiplicative 
cascade the smaller k, and conversely [Fig. 2b in [52]]. 

One can explore the Tsallis type of the probability distribution function Eq.(4) in two 
limits. For small values of normalized log variability Z the probability distribution function 
converges to the form 



small values of the normalized log variability, for any At investigated hereby (sec Figs. 2-3). 

In the limit of large values of normalized variability Z, the Tsallis type distribution 
converges to a power law 



Studying the Tsallis type of distribution function one can obtain from Eq.(4) an ex- 
pression for the width of the TsaUis type of probabihty distribution function, 2(7^ = 
{2a - {q - l))/{2aCPo{q - !))■ In the limit of a — > 1 the width of the Tsallis type 
distribution 2cr^ = (3 — q)/2C/3o{q — 1), i.e. ~ 2/(C/3o)- It is obvious that for large time 
lags 2(7^ tends to diverge [32], like ~ (At)^/(^~^); this can be easily verified on a log-log plot 
(not shown). 




(14) 



Therefore the Tsallis type distribution function converges to a Gaussian, i.e. a 



1, for 




(15) 
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In limit of q — > 1 the Tsallis type distribution function converges to Gaussian. The 
values of the parameters g, a, C/3o, that best fit the data using Eq.(4), and 2(7^ are plotted 
as a function of the time lag in Fig. 10. 

V. FOKKER-PLANCK APPROACH 

On the other hand, the evolution of a time dependent probability distribution function is 
usually described within the Fokker-Planck approach. This method provides some further 
information on the correlations present in the time series and it begins with the joint PDF's, 
that depend on jV variables, i.e. p^{Zi, Ati; Z/j-, At^). We started to address this issue 
by determining the joint PDF for J\f = 2, i.e. p(Z2, At2; Axi, Ati). The symmetrically 
tilted character of the joint PDF contour levels (Fig. 11) around an inertia axis with slope 
+1 points out to some statistical dependence, i.e. a correlation, between the normalized 
variabihty Z{t, At) of the monthly value signal of SOI. 

The conditional probability function is 

piZ.,MZ.,M,>^ ^('-^'^;^'^-^ (16, 

for i — 1, ...,jV — 1. For any At2 < Atj < A^i, the Chapman-Kolmogorov equation is a 
necessary condition of a Markov process, one without memory but governed by probabilistic 
conditions 

p{Z2, At2\Z,, Ah) = J d{Zi)p{Z2, At2\Zi, Ati)p{Axi, AU\Zi, Ah). (17) 

The Chapman-Kolmogorov equation when formulated in differential form yields a mas- 
ter equation, which can take the form of a Fokker-Planck equation [43]. Let r = log2{^S/ At), 



i2 



p{Z,T) (18) 



in terms of a drift D^^\Z,t) and a diffusion coefficient D^'^\Z,t) (thus values of r represent 
Ati, i = I,-)- 

The coefficient functional dependence can be estimated directly from the moments M*^^) 
(known as Kramers-Moyal coefficients) of the conditional probability distributions: 
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M^'^) = ^ y dZ'{Z' - Zfv{Z\ T + At|Z, t) (19) 

D('=)(Z,t) = ^limM^ (20) 

for At 0. The drift coefficient D^^^ and the diffusion coefficient are well represented 
(Fig. 12a,b) by a line and parabola, respectively 

= -0.37Z - 0.01 (21) 

= 0.10^2-0.10^ + 0.33 (22) 

for the normalized variability (plotted with dots). 

We have compared the above values of the drift and the diffusion coefficients for those 
of the drift and the diffusion coefficients for non normalized variability (plotted with open 
circles) and have obtained 

= -0.52Ay - 0.02 (23) 

£1(2) = 0.20A2/2 - O.lOAy + 0.24 (24) 

Note that the first term on the right hand side of Eq. (18) is identified [54] as the term 
generating drift behavior in the evolution of the PDF, while the second term is responsible 
for the diffusion, or fluctuation term in the PDF evolution. In the asymptotic case when 
the linear term in D*^^), i.e. the coefficient D^^^ is dominating the dependence and the 
independent Dq^'* term is somewhat dominating in D^^^ then the Fokker-Planck equation 
is linear, otherwise the drift and diffusion terms are intervened. Comparing the values 
of the linear and independent terms in Eqs. (21) and (23), as well as the values of the 
quadratic, linear and independent terms in Eqs. (22) and (24) one may argue that the 
linear approximation for the Fokker-Planck equation holds more convincing for the PDF 
evolution of the normalized variability (Eqs. (21,23)) as opposite to the PDF evolution of 
the non normalized variability Ay. 
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On the other hand, it may be worthwhile to recall that the observed quadratic depen- 
dence of the diffusion term D^^^ is essential for the logarithmic scahng of the intermittency 
parameter in studies on turbulence. 

Finally, the Fokker-Planck equation for the distribution function is known to be equivalent 
to a Langevin equation for the variable, i.e. Z here, (within the Ito interpretation [43-46, 53]) 

^Z(r) = i?«(Z(r),r)+77(r)^L>(2)(Z(r),r), (25) 

where r]{T) is a fluctuating ^-correlated force with Gaussian statistics, i.e. < ?7(r) ri{T')> = 
2S{t-t'). 

Thus the Fokker-Planck approach provides the evolution process of PDF's from small 
time lags to larger ones. The fact that the drift coefficient is finite implies that there is 
some 'restoring force', i.e. 7 7^ in Eq. (6), while the quadratic dependence of D^^^ in Z is 
obviously like an autocorrelation function for a diffusion process. 

An interaction that can produce such a 'restoring force' is the the air-sea interaction that 
takes place in many different ways. However, its main ingredients are the air-sea fluxes of 
mass, heat and momentum. Reliable estimates of the air-sea fluxes of heat and momentum 
are vital to improve our understanding of the coupled ocean-atmosphere system. An air-sea 
heat and momentum climatology, as the one recently reported [55] can be hopefully used 
in future dynamic models to relate the flndings of this study to classical meteorological 
qualities using observations. 

VI. CONCLUSION 

In summary, we have presented a method that provides the evolution process of proba- 
bility distribution functions (over 140 years) for one climatological index, i.e. the SOI. It 
can be flrst recalled that crossovers in the DFA results and power spectral density point to 
speciflc time scales, - in fact related to famous phenomena, like sunspots. 

We have mainly studied the evolution process of the tails that are outside the central 
(Gaussian) regime,- which of course occurs only at small variability, thereby facilitating the 
understanding of the evolution of these distribution functions in a Fokker-Planck framework. 
The Gaussian regime range has been found to imply that signal correlations extend up to 
ca. 48 months, - an interesting time lag to be considered in microscopic evolution model(s) 
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of El Nino. 

It has been found that Beck turbulence model can be well applied to describe the dis- 
tributions of the standard deviation of the SOI signal normalized variability assuming a 

X^-distribution for these. In some sense this application (or generalization of the model) 
is justified by the fact that the ideas behind the turbulence model, based on temperature 
fluctuations, can be expected to be carried over to the case in which pressure fluctuations 
occur in the system. 

An open question in nonextensive thermostatistics studies is often raised about the mean- 
ing, value and behavior of the non extensive exponent, or Tsallis parameter q. The inter- 
mittency exponent is interestingly found to be related to the scaling exponent of the PDF 
moments in the framework of Kolmogorov log-normal model, thereby giving weight to the 
model and the statistical approach. 

We have also presented the turbulence-like dynamics through the Fokker- Planck and the 
Langevin equations. We have (as it has been expected) found that, in the treated case, 
there is some 'restoring force', i.e. (7 7^ in the Langevin equation). A comparison is made 
between normalized variability and non normalized variability. 

Whence we have related a climatological signal behavior to Tsallis non extensive ther- 
modynamics approach, i.e. more precisely to a turbulence-like process, - as climatological 
ocean- atmospheric interface interactions and indices were often claimed to be seen. No need 
to say that this empirical modeling only describes the evolution of the signal but does not 
explain it, as a general circulation model [56] should do. Nevertheless the time scales which 
are hereby observed might shine some light on approximation validity or the need to restrict 
extrapolations to realistic ranges, - including memory effects. Finally, it seems that we have 
thoroughly answered the often raised question 'why to look at the tails of a probability 
distribution function? and what does that lead to?'. 



Acknowledgements 

Part of FP work has been supported by European Commission Project E2C2 FP6-2003- 

NEST-Path-012975 Extreme Events: Causes and Consequences. Part of this work results 
from flnancing through the ARC 02-07/293 Project of the ULg and the COST PIO "Physics 



15 



of Risk" program which MA also thoroughly acknowledges. Critical and encouraging com- 
ments by M. Ghil, H. Herrmann, K. Ivanova, J. Peinke, C. TsaUis, and J. Vannimenus have 
been as always very valuable for improving this report. 



[1] B.B. Mandelbrot and J.W. Van Ness, SIAM Rev., 10, 422 (1968) 

[2] K. Ivanova, M. Ausloos, E.E. Clothiaux and T.P. Ackerman, Europhys. Lett. 52, 40 (2000) 
[3] E. Kosciclny-Bunde, A. Bunde, S. Havlin, H. E. Roman, Y. Goldreich and H.-J. Schellnhuber, 

Phys. Rev. Lett. 81, 729 (1998) 
[4] D.L. Turcotte, Fractals and Chaos in Geology and Geophysics (Cambridge University Press, 

Cambridge, 1997) 

[5] P. Ch. Ivanov, M. G. Rosenblum, C.-K. Peng, J. E. Mietus, S. Havlin, H. E.Stanley and 

A. L. Goldbergcr, Nature 383, 323 (1996) 
[6] K. Ivanova and M. Ausloos, Eur. Phys. J. B 27, 239 (2002) 
[7] N. Vandewalle and M. Ausloos, Physica A 246, 454 (1997) 
[8] N. Vandewalle and M. Ausloos, Phy. Rev. E 58, 6832 (1998) 
[9] P. Flandrin, IEEE Trans. Inform. Theory 35, 197 (1989) 
[10] F. M. Ramos, C. Rodrigucs Ncto, R. R. Rosa, L. D. Abreu Sa and M. J. A. Bolzan, Nonlinear 

Analysis - Theory 23, 3521 (2001) 
[11] U. Prisch, Turbulence: the Legacy of A. N. Kolmogorov, (Cambridge Univ. Press, Cambridge, 

1995) 

[12] C. Tsallis, J. Stat. Phys. 52, 479 (1988) 

[13] C. Tsallis and D.J. Biikman, Phys. Rev. E 54, 2197 (1996) 

[14] G. Wilk and Z. Wlodarczyk, Phys. Rev. Lett. 84, 2770 (2000) 

[15] T. Arimitsu and N. Arimitsu, Prog. Theor. Phys. 105, 355 (2001) 

[16] C. Beck, Phys. Rev. Lett. 87, 180601 (2001) 

[17] C. Beck, Physica A 295, 195 (2001) 

[18] C. Beck, G. S. Lewis, and H. L. Swinney, Phys. Rev. E 63, 035303 (2001) 



16 



[19] F. Sattin, Physica Scripta 71, 443 (2005) 
[20] F. Sattin, Phys. Rev. E 68, 032102 (2003) 

[21] A report on the matter was presented by the authors at the 2006 EGU General Assembly, 

Vienna, April 06. 
[22] O.K. Vallis, J. Geophys. Res. 93 13 979 (1988). 
[23] S.G. Philander, Rep. Prog. Phys. 62, 123 (1999). 
[24] M. Ghil and N. Jiang Geophys. Res. Lett. 25, 171 (1998). 
[25] C.L. Kcppenne and M. Ghil J. Geophys. Res. 97, 20449 (1992). 
[26] M. Ausloos and K. Ivanova , Phys. Rev. E 63, 047201 (2001). 

[27] R. Garcia, P. Ribera, L. Gimenoo, E. Hernandez Ann. Geophysicae 18, 247 ( 2000). 
[28] C. Collette and M. Ausloos, Int. J. Mod. Phys. C 15, 1353 (2005) 
[29] G.T. Walker, Q. J. R. Meteorol. Soc. 54, 79 (1928). 

[30] P.J. Brockwell and R.A. Davis, Time Series : Theory and Methods (Springer- Verlag, Berlin, 
1991) 

[31] http : // www.hom.gov.au/ climate/ glossary / soi.shtml 
[32] F. Michael and M.D. Johnson, Physica A 320, 525 (2003). 
[33] N Kozuki and N Fuchikami, Physica A 329, 222 (2003). 

[34] G.P. Konncn, P.D. Jones, M.H. Kaltofcn and R.J. Allan, J. Climate, 11, 2325 (1998) 
[35] C. Tsallis, C. Anteneodo, L. Borland and R. Osorio, Physica A 324, 89 (2003). 
[36] M.S. Taqqu, V. Teverovsky and W. Willinger, Fractals 3, 785 (1995) 

[37] K. Hu, Z. Chen, P. Ch. Ivanov, P. Carpena and H.E. Stanley, Phys. Rev E 64, 011114 (2001) 
[38] K. Ivanova and M. Ausloos, Physica A 274, 349 (1999) 

[39] G.M. Raymond, D.B. Percival and J.B. Bassingthwaighte, Physica A 322, 169 (2003) 
[40] T.C. Halsey, M.H. Jensen, L.P. Kadanoff, I. Procaccia and B.I. Shraiman, Phys. Rev. A 33, 
1141 (1986) 

[41] C. Tsallis, S.V.F. Levy, A.M.C. Souza and R. Maynard, Phys. Rev. Lett. 75, 3589 (1995) 
[42] F. Sattin, J. Phys. A 36, 1583 (2003) 

[43] M.H. Ernst, L.K. Haines and J.R. Dorfman, Rev. Mod. Phys. 41, 296 (1969) 

[44] L.E. Reichl, A Modern Course in Statistical Physics, Univ. Texas Press, Austin (1980) 

[45] P. Hanggi and H. Thomas, Phys. Rep. 88, 207 (1982) 

[46] C.W. Gardiner, Handbook of Stochastic Methods, Springer- Verlag, Berlin (1983) 

17 



[47] A.M. Kolmogorov, J. Fluid Mech. 13, 82 (1962) 

[48] K.R. Sreenivasan and P. Kailasnath, Phys. Fluids 5, 512 (1993) 

[49] K.R. Sreenivasan, R.A. Antonia and H.Q. Danh, Phys. Fluids 20, 1238 (1977) 

[50] R.R. Prasad, C. Meneveau and K.R. Sreenivasan, Phys. Rev. Lett 61, 74 (1988) 

[51] J.C. Wyngaard and H. Tennekes, Phys. Fluids 13, 273 (1970) 

[52] R.A. Janik and B. Ziaja, Acta Phys. Pol. B 30, 259 (1999) 

[53] H. Risken, The Fokker-Planck Equation: Methods of Solution and Applications, 2nd edn, 
(Springer- Verlag, Berlin, 1989) 

[54] H.J. Carmichael, Statistical Methods in Quantum Optics 1. Master Equations and Fokker- 
Planck Equations, (Springer- Verlag, Berlin, 1999). 

[55] S.A. Josey, E.G. Kent and P.K.Taylor, Southampton Oceanography Centre Report, 6, 30 (1999) 

[56] http : / /www — das.uwyo.edu/ geerts / cwx / notes / chapl2/ nwpgcm.html 

[57] Sattin formula [20] might also be tested in future work. 



18 



Figure Captions 



Figure 1 - Monthly values of the Southern Oscillation Index (SOI) as defined in 
the text reported from Jan. 1866 to Jan. 2006. Data are downloaded from 
http : / /tao.atmos.washington.edu/pacs/additionalanalyses/soi.html and from ftp : 
// ftpprd.ncep.noaa.gov/pub/cpc/wd52dg/data/indices/soi after June 1999. Data series 
consists of 1681 data points 

Figure 2 - (a) Probability distribution function of normalized variability Z{t, At) of 
monthly values signal of the Southern Oscillation Index from Jan. 1866 to Jan. 2006 
for At ^ 1 month (symbols). Z{t,At) is defined as Z(t,At) = {Ay{t)- < Ay >At)/<^At, 
where Ay{t) — y{t + At) — y{t) and aAt is the standard deviation of Ay{t) for time lag At. 
The dashed fine represents a Gaussian distribution. Inset: Power law fit (solid fine) of the 
negative and positive slope (-5.6 for both) of the distribution outside the Gaussian regime, 
i.e. ] — 2, +2[. (b) same as (a) but for non normalized variability Ay 

Figure 3 - Probability distribution function (PDF) pAt{Z) of normalized variability of 
monthly values signal of the SOI (symbols) and the Tsallis type distribution function (lines) 
for different values of At — 1, 3, 6, 12, 24, 36 months. The PDF (symbols and curves) for each 
At are moved down by 10 with respect to the previous one; the curve for At — 1 month is 
unmoved. The large dots mark the ends of the interval in which the distribution is like a 
gaussian distribution. The values of the parameters for the Tsallis type distribution function 
for each At are summarized in Table 11 

Figure 4 - DFA function F{n) plotted as a function of the the box size n of the integrated 
monthly values signal of the SOI from Jan. 1866 to Jan. 2006. l//-like fluctuations with 
slope = 1.06 ± 0.01 are obtained for time scales below 66 months and fractional Gaussian 
noise like fluctuations slope — 0.36 ± 0.02 above 72 months. Insets: White noise like 
fluctuations of two types of surrogate data, when the data are shuffled randomly and when 
the sign of the data is shuffled randomly 

Figure 5 - DFA function F{n) plotted as a function of the box size n of the integrated 
normalized variability Z{t, At) of the monthly values signal of the SOI from Jan. 1866 
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to Jan. 2006, for different time lags At — 1,3,6,12,24,36 montfis. Values of the scaling 
exponents Hdfa for the various DFA functions are summarized in Table I 

Figure 6 - Power spectrum S{f) of the monthly values signal of the SOI from Jan. 1866 to 
Jan. 2006. A scale break at around / = 1/70 month"^ separates two scahng regions. Insets: 
Scaling of the power spectrum of both shuffled amplitude and shuffled sign of monthly values 
signal of the SOI as a white noise signal with a; 

Figure 7 - Power spectrum S{f ) of the normalized variability Z(t, At) of the monthly 
values signal of the SOI from Jan. 1866 to Apr. 2003 for different time lags At = 
1, 3, 6, 12, 24, 36 months. Each curve is moved down by 10~^ with respect to the previous 
one; the power spectrum of the normalized returns for At — 1 month is not displaced 

Figure 8 - Probability density fAt{P) of the local volatihty P (Eq.(lO)) in terms of standard 
deviation of the normalized variability Z{t, At) of SOI in non-overlapping windows with size 
m=12 months for different time lags (symbols) (a-f ) At — 1, 3, 6, 12, 24, 36 months. Lines: 
X^-distribution as given by Eq. (7) 

Figure 9 - The functional dependence of the Tsallis q parameter on the rescaled time lag 
At/Ati for AtL = 48 months and 6 = 0.11 (see Eq. (13)) (line); the symbols represent 
the values of the q parameter listed in Table II and used to plot the fitting lines in Fig. 2. 
Inset : Scahng properties of the rescaled kurtosis K^/Kl, where = 3 is the kurtosis for 
a Gaussian process, as a function of the rescaled time lag At/ Ati satisfying Eq. (11) (open 
symbols) and Eq. (12) (full symbols) 

Figure 10 - Characteristic parameters of TsaUis type distribution function as defined in 
[33] : Tsallis g-parameter (crosses), a (squares), constant C/3o used in the fit (open circles), 
the width of the Tsallis type distribution 2(T^ = {2a — {q — 1)) / {2aC j3Q{q — 1)) from Eq.(4) 
(triangles) (rescaled by a factor of 1/180, for better display) 

Figure 11 - Typical contour plots of the joint probability density function 
p{Z2, At2; Zi, Ati) of the monthly values signal of the SOI from Jan. 1866 to Jan. 2006. 
Dashed lines have a slope +1 and emphasize the correlations between probability den- 
sity functions for At2 — 1 month and Ati — 2 months. Contour levels correspond to 



20 



logiop{Z2, Aia; Zi, Ati) = -1.8, -2.0, -2.2, -2.4, -2.6, -2.8 from center to border 

Figure 12 - Kramers- Moyal drift (a) D^^^ and diffusion (b) D^'^^ coefficients as a function 
of normalized variability Z (dots) and non normalized variability Z (open circles) of the 
monthly values signal of the SOI; L>W = -0.37Z-0.01 ((a) dots), D^^) = 0.10^^-0.10^ + 
0.33 ((b) dots); L>W = -0.52A|/ - 0.02 ((a) open circles), D^^^ = 0.20A|/2 - O.IOA7/ + 0.24 
((b) open circles) 
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TABLE I: Values of the scaling exponent from the DFA analysis of normalized variability Z for 
different values of the time lag Ai = 1, 3, 6, 12, 24, 36 months, and crossover 'box size' Ux 

At 1 + Hdfai 1 + HdfA2 

1 0.264± 0.017 0.058± 0.016 23 

3 0.594± 0.038 0.082± 0.060 23 

6 0.909± 0.028 O.lOli 0.063 23 
12 1.148± 0.032 0.129± 0.070 25 
24 1.115± 0.030 0.142± 0.066 32 
36 1.082± 0.027 0.192± 0.044 35 



TABLE II: Values of the parameters characterizing the monthly values of the Southern Oscillation 
Index (SOI) from Jan. 1866 to Jan. 2006 in the nonextensive thermostatistics approach. For the 
definition of the Kolmogorov-Smirnov distance (Iks see the text 
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FIG. 1: Monthly values of the Southern Oscillation Index (SOI) as defined in the 
text reported from Jan. 1866 to Jan. 2006. Data are downloaded from 
http : / /taoMtmos.washington.edu/pacs/additionalanalyses/ soi.html and from ftp : 
// f tpprd.ncep.noaa.gov/pub/cpc/wd52dg/data/indices/ soi after June 1999. Data series consists 
of 1681 data points 
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non normalized variability 



FIG. 2: (a) Probability distribution function of normalized variability Z{t, At) of monthly values 
signal of the Southern Oscillation Index from Jan. 1866 to Jan. 2006 for At = 1 month (symbols). 
Z{t,At) is defined as Z{t,At) = {Ay{t)- < Ay >At)/crAt, where Ay{t) = y{t + At) - y{t) and 
aAt is the standard deviation of Ay(t) for time lag At. The dashed line represents a Gaussian 
distribution. Inset: Power law fit (solid line) of the negative and positive slope (-5.6 for both) of 
the distribution outside the Gaussian regime, i.e. ] — 2, +2[. (b) same as (a) but for non normalized 
variability Ay 
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normalized variability Z 



FIG. 3: Probability distribution function (PDF) pAt{Z) of normalized variability of monthly values 
signal of the SOI (symbols) and the Tsallis type distribution function (lines) for different values of 
At = 1, 3, 6, 12, 24, 36 months. The PDF (symbols and curves) for each At are moved down by 10 
with respect to the previous one; the curve for At = 1 month is unmoved. The large dots mark 
the ends of the interval in which the distribution is like a gaussian distribution. The values of the 
parameters for the Tsallis type distribution function for each At are summarized in Table II 
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FIG. 4: DFA function F{n) plotted as a function of the the box size n of the integrated monthly 
values signal of the SOI from Jan. 1866 to Jan. 2006. l//-like fluctuations with slope = 1.06ib0.01 
are obtained for time scales below 66 months and fractional Gaussian noise like fluctuations slope = 
0.36 lb 0.02 above 72 months. Insets: White noise like fluctuations of two types of surrogate data, 
when the data are shuffled randomly and when the sign of the data is shuffled randomly 
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FIG. 5: DFA function F{n) plotted as a function of the box size n of the integrated normahzed 
variabihty Z{t,At) of the monthly values signal of the SOI from Jan. 1866 to Jan. 2006, for 
different time lags At = 1,3,6,12,24,36 months. Values of the scaling exponents Hdfa for the 
various DFA functions are summarized in Table I 
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FIG. 6: Power spectrum S{f) of the monthly values signal of the SOI from Jan. 1866 to Jan. 2006. 
A scale break at around / = 1/70 month~^ separates two scaling regions. Insets: Scaling of the 
power spectrum of both shuffled amplitude and shuffled sign of monthly values signal of the SOI 
as a white noise signal with ^ ~ 
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FIG. 7: Power spectrum S{ f) of the normalized variability Z{t, At) of the monthly values signal 
of the SOI from Jan. 1866 to Apr. 2003 for different time lags At = 1,3,6,12,24,36 months. 
Each curve is moved down by 10~^ with respect to the previous one; the power spectrum of the 
normalized returns for At = 1 month is not displaced 
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FIG. 8: Probability density fAt{P) of the local volatility /? (Eq.(lO)) in terms of standard deviation 
of the normalized variability Z{t, At) of SOI in non-overlapping windows with size m=12 months 
for different time lags (symbols) (a-f) At = 1,3,6,12,24,36 months. Lines: x^"distribution as 
given by Eq. (7) 
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+ fit in Fig. 2 




FIG. 9: The functional dependence of the Tsallis q parameter on the rescaled time lag At/AtL for 
= 48 months and S = 0.11 (see Eq. (13)) (line); the symbols represent the values of the q 
parameter listed in Table II and used to plot the fitting lines in Fig. 2. Inset : Scaling properties 
of the rescaled kurtosis Kr/Ki, where Kl = 3 is the kurtosis for a Gaussian process, as a function 
of the rescaled time lag At/Ati satisfying Eq. (11) (open symbols) and Eq. (12) (full symbols) 
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FIG. 10: Characteristic parameters of Tsallis type distribution function as defined in [33] : Tsallis 
g-parameter (crosses), a (squares), constant C/?o used in the fit (open circles), the width of the 
Tsalhs type distribution 2o-^ = (2a — {q — l))/(2aC/3o(^ — 1)) from Eq.(4) (triangles) (rescaled by 
a factor of 1/180, for better display) 
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FIG. 11: Typical contour plots of the joint probability density function p{Z2, At2; Zi, Ati) of 
the monthly values signal of the SOI from Jan. 1866 to Jan. 2006. Dashed lines have 
a slope +1 and emphasize the correlations between probability density functions for At2 = 
1 month and Ati = 2 months. Contour levels correspond to logiop{Z2, At2; Zi, Ati) = 
-1.8, -2.0, -2.2, -2.4, -2.6, -2.8 from center to border 
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FIG. 12: Kramers- Moyal drift (a) D^^^ and diffusion (b) D*-^-* coefficients as a function of normalized 
variability Z (dots) and non normalized variability Z (open circles) of the monthly values signal 
of the SOI; L>(i) = -0.37Z - 0.01 ((a) dots), D^^) ^ o.lOZ^ - O.IOZ + 0.33 ((b) dots); D^^^ = 
-0.52A?/ - 0.02 ((a) open circles), D^^) = o.20Ay2 _ O.lOAy + 0.24 ((b) open circles) 
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